We include some more details about the \(\texttt{Diel.Niche}\) package. For the fundamental uses, see the main vignette (’Diel-Niche-vignette).
Along with the three main hypotheses sets (Maximizing, Traditioanl, and General), we have two additional sets. These sets are importantly different in that the entire parameter space is not defined. They are also defined a bit differently.
Threshold
This hypothesis set is based on defining each diel period in terms of lower thresholds. If the thresholds are not provided, default values are used.
plot.diel(hyp=hyp.sets("Threshold"))Not that where there is no colored space indicates that this set of probabilities are not considered as possible parameter values.
Variation
This hypothesis set is based on defining the ranges (lower, upper) of possible probability values.
plot.diel(hyp=hyp.sets("Variation"))You do not need to rely on pre-set definitions of diel hypotheses. If you are interested in categorizing dominant diel phenotypes based on 0.90 probability or higher, you can adjust the hypotheses using the function ‘diel.ineq()’ than pass this object to the plot function or model fitting function. Here, we will just define the lower bount of each diurnal, nocturnal, and crepuscular using \(\xi_{1} = 0.90\)
diel.setup=diel.ineq(xi = c(0.90))
plot.diel(hyp=hyp.sets("Traditional"),diel.setup=diel.setup)If you want to separate the hypotheses and reduce the parameter space of cathemerality,
diel.setup=diel.ineq(xi = 0.70,separation=0.10)
plot.diel(hyp=hyp.sets("Traditional"),diel.setup=diel.setup)When considering the General hypothesis set, you can adjust the lower probability value (\(\xi_{1}\)) to define crepuscular, diurnal, and nocturnal, as well as adjust the definitions of Cathemeral General vs Diurnal-Nocturnal, Crepuscular-Nocturnal, and Diurnal-Crepuscular (\(\xi_{2}\)). To do this, set \(xi\) as a vector of two probabilities,
diel.setup=diel.ineq(xi = c(0.9,0.05))
plot.diel(hyp=hyp.sets("General"),diel.setup=diel.setup)This new hypothesis set can be used to simulate data and fit the same models by passing the object diel.setup as an attribute of the functions sim.diel and diel.fit.
y=sim.diel(hyp="D.CR",n.sample = 200,diel.setup=diel.setup)$y
out = diel.fit(y,hyp.set=hyp.sets("General"),n.chains=2,post.fit = TRUE,diel.setup=diel.setup)
#> Data checks Complete.
#> Calculating Bayes Factors...
#> Posterior Sampling...
#> The most supported model is:
#> Cathemeral (General)
plot.diel(out)There are additional defaults that could be changed, depending on your hypotheses of interest.
#To see all the default inputs
defaults=diel.ineq()
defaults$inputs
#> $e.D
#> [1] 0.1
#>
#> $e.N
#> [1] 0.1
#>
#> $e.CR
#> [1] 0.1
#>
#> $e.EC
#> [1] 0.1
#>
#> $e.AV
#> [1] 0.1
#>
#> $xi.D
#> [1] 0.8 0.9
#>
#> $xi.N
#> [1] 0.8 0.9
#>
#> $xi.CR
#> [1] 0.8 0.9
#>
#> $xi.C
#> [1] 0.2
#>
#> $xi
#> [1] 0.8 0.1
#>
#> $xi.EC
#> [1] 0.33
#>
#> $p.avail
#> [1] 0.1666660 0.4166667
#>
#> $separation
#> [1] 0For the Variation hypotheses, we can adjust the range of values for all hypotheses as,
# Default
diel.plot(hyp =hyp.sets("Variation"))For the Threshold hypotheses, we can adjust the lower probability similarly as,
# Default
diel.plot(hyp =hyp.sets("Threshold"))Next, we will demonstrate how to include additional diel hypotheses. This requires understanding the inequality constraint setup and then coding the results as additions to the ‘diel.ineq’ output. Let’s consider two complementary hypotheses that account for the complete parameter space among the three fundamental probabilities. Our objective is to evaluate whether an animal is highly crepuscular or not. We will define highly crepuscular as \(p_{tw} \geq 0.5\). The complementing hypothesis is that the animal does not have such high crepuscular activity \(p_{tw} < 0.5\).
\[ \begin{equation} \begin{split} & p_{\text{tw}} \geq 0.5 \\ \end{split} \end{equation} \]
Next, we need to translate each inequality into the framework of \(\mathbf{A} \boldsymbol{\theta} \leq \mathbf{b}\).
\[ \begin{equation} \begin{split} & p_{\text{tw}} \geq 0.50 \\ & -p_{\text{tw}} \leq -0.50 \\ & (-1) \times p_{\text{tw}} + (0)\times p_{\text{d}} \leq -0.50 \\ \end{split} \end{equation} \] Putting it together, we take the constants in parentheses left of the equal sign and package it into matrix \(\mathbf{A}\) and take the constants on the right of the equal sign and package it into vector \(\mathbf{b}\) as,
\[ \begin{equation} \textbf{A} = \begin{bmatrix} -1 & 0 \\ \end{bmatrix}, \textbf{b} = \begin{bmatrix} -0.50 \\ \end{bmatrix}. \end{equation} \]
The complementing low activity during twilight as,
\[ \begin{equation} \begin{split} & p_{\text{tw}} \leq 0.49999 \\ \end{split} \end{equation} \] Following the same procedure, we find that,
\[ \begin{equation} \textbf{A} = \begin{bmatrix} 1 & 0 \\ \end{bmatrix}, \textbf{b} = \begin{bmatrix} 0.49999 \end{bmatrix}. \end{equation} \] ## Coding
Now, we need to code these hypotheses into \(\texttt{Diel.Niche}\)
# First, we need to create new list elements. One for each new hypothesis, with these ordered elements in each list,
# Element 1: Name
# Element 2: A
# Element 3: b
# Element 4: func
# Note, func - this should almost always be 'bf_multinom', unless an equality statement is specified, which required function 'bf_equality'
A <- matrix(c(-1,0),ncol = 2, byrow = TRUE)
b <- c(-0.5)
New1=list(Name="Highly Crepuscular",A=A,b=b,func="bf_multinom")
A <- matrix(c(1,0),ncol = 2, byrow = TRUE)
b <- c(0.49999)
New2=list(Name="Not highly Crepuscular",A=A,b=b,func="bf_multinom")
#Second, include these new lists into diel.ineq() function as,
diel.setup=diel.ineq()
diel.setup$New1=New1
diel.setup$New2=New2We have setup our hypotheses. Now, we need to visualize them to check it’s what we intended. Note that in each following, we pass our new diel.setup object to each function.
plot.diel(hyp=c("New1","New2"),diel.setup = diel.setup,more.points = TRUE)Let’s simulate data under the hypothesis that the animal is highly crepuscular.
set.seed(43243)
dat=sim.diel(n.sample=85, hyp="New1",diel.setup = diel.setup)
y=dat$y
y
#> y_crep y_day y_night
#> [1,] 61 22 2Next, let’s compare our two models
#Note, that we get a warning, which is fine
out <- diel.fit(y, hyp.set = c("New1","New2"),diel.setup=diel.setup)
#> Data checks Complete.
#> Calculating Bayes Factors...
#> Warning in print.hyp.name(ms.model): Hypothesis code is not recognized
#> The most supported model is:
#> New1
#Let's look at the model probabilities, where we find the model we simulated is highly supported
out$bf.table
#> Prior Posterior
#> New1 0.5 1
#> New2 0.5 0